home *** CD-ROM | disk | FTP | other *** search
/ Die Ultimative Software-P…i Collection 1996 & 1997 / Die Ultimative Software-Pakete CD-ROM fur Atari Collection 1996 & 1997.iso / g / gnu_c / pmlsrc23.zoo / pmlsrc / asin.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-03-19  |  5.3 KB  |  269 lines
  1. /************************************************************************
  2.  *                                    *
  3.  *                N O T I C E                *
  4.  *                                    *
  5.  *            Copyright Abandoned, 1987, Fred Fish        *
  6.  *                                    *
  7.  *    This previously copyrighted work has been placed into the    *
  8.  *    public domain by the author (Fred Fish) and may be freely used    *
  9.  *    for any purpose, private or commercial.  I would appreciate    *
  10.  *    it, as a courtesy, if this notice is left in all copies and    *
  11.  *    derivative works.  Thank you, and enjoy...            *
  12.  *                                    *
  13.  *    The author makes no warranty of any kind with respect to this    *
  14.  *    product and explicitly disclaims any implied warranties of    *
  15.  *    merchantability or fitness for any particular purpose.        *
  16.  *                                    *
  17.  ************************************************************************
  18.  */
  19.  
  20.  
  21. /*
  22.  *  FUNCTION
  23.  *
  24.  *    asin   double precision arc sine
  25.  *
  26.  *  KEY WORDS
  27.  *
  28.  *    asin
  29.  *    machine independent routines
  30.  *    trigonometric functions
  31.  *    math libraries
  32.  *
  33.  *  DESCRIPTION
  34.  *
  35.  *    Returns double precision arc sine of double precision
  36.  *    floating point argument.
  37.  *
  38.  *    If argument is less than -1.0 or greater than +1.0, calls
  39.  *    matherr with a DOMAIN error.  If matherr does not handle
  40.  *    the error then prints error message and returns 0.
  41.  *
  42.  *  USAGE
  43.  *
  44.  *    double asin (x)
  45.  *    double x;
  46.  *
  47.  *  REFERENCES
  48.  *
  49.  *    Fortran IV-plus user's guide, Digital Equipment Corp. pp B-2.
  50.  *
  51.  *  RESTRICTIONS
  52.  *
  53.  *    For precision information refer to documentation of the floating
  54.  *    point library primatives called.
  55.  *    
  56.  *  PROGRAMMER
  57.  *
  58.  *    Fred Fish
  59.  *
  60.  *  INTERNALS
  61.  *
  62.  *    Computes arcsine(x) from:
  63.  *
  64.  *        (1)    If x < -1.0 or x > +1.0 then
  65.  *            call matherr and return 0.0 by default.
  66.  *
  67.  *        (2)    If x = 0.0 then asin(x) = 0.0
  68.  *
  69.  *        (3)    If x = 1.0 then asin(x) = PI/2.
  70.  *
  71.  *        (4)    If x = -1.0 then asin(x) = -PI/2
  72.  *
  73.  *        (5)    If -1.0 < x < 1.0 then
  74.  *            asin(x) = atan(y)
  75.  *            y = x / sqrt[1-(x**2)]
  76.  *
  77.  */
  78.  
  79. #include <stdio.h>
  80. #include <math.h>
  81. #include "pml.h"
  82.  
  83. #if !defined (__M68881__) && !defined (sfp004)    /* mjr++    */
  84.  
  85. static char funcname[] = "asin";
  86.  
  87.     struct exception xcpt;
  88.  
  89. double asin (x)
  90. double x;
  91. {
  92.  
  93.     if ( x > 1.0 || x < -1.0) {
  94.     xcpt.type = DOMAIN;
  95.     xcpt.name = funcname;
  96.     xcpt.arg1 = x;
  97.     if (!matherr (&xcpt)) {
  98.         fprintf (stderr, "%s: DOMAIN error\n", funcname);
  99.         errno = EDOM;
  100.         xcpt.retval = HUGE_VAL;
  101.     }
  102.     } else if (x == 0.0) {
  103.     xcpt.retval = 0.0;
  104.     } else if (x == 1.0) {
  105.     xcpt.retval = HALFPI;
  106.     } else if (x == -1.0) {
  107.     xcpt.retval = -HALFPI;
  108.     } else {
  109.     xcpt.retval = atan ( x / sqrt (1.0 - (x * x)) );
  110.     }
  111.     return (xcpt.retval);
  112. }
  113. #endif !defined (__M68881__) && !defined (sfp004)    /* mjr++    */
  114. #ifdef    sfp004
  115.  
  116. __asm("
  117.  
  118. comm =     -6
  119. resp =    -16
  120. zahl =      0
  121.  
  122. ");    /* end asm    */
  123.  
  124. #endif    sfp004
  125. #if defined (__M68881__) || defined (sfp004)
  126.  
  127.     __asm(".text; .even");
  128.  
  129. # ifdef ERROR_CHECK
  130.  
  131.     __asm("
  132.  
  133. _Overflow:
  134.     .ascii \"OVERFLOW\\0\"
  135. _Domain:
  136.     .ascii \"DOMAIN\\0\"
  137. _Error_String:
  138.     .ascii \"asin: %s error\\n\\0\"
  139. .even
  140.  
  141. | pml compatible asingent
  142. | m.ritzert 7.12.1991
  143. | ritzert@dfg.dbp.de
  144. |
  145. |    /* NAN  = {7fffffff,ffffffff}        */
  146. |    /* +Inf = {7ff00000,00000000}        */
  147. |    /* -Inf = {fff00000,00000000}        */
  148. |    /* MAX_D= {7fee42d1,30773b76}        */
  149. |    /* MIN_D= {ffee42d1,30773b76}        */
  150.  
  151. .even
  152. double_max:
  153.     .long    0x7fee42d1
  154.     .long    0x30273b76
  155. double_min:
  156.     .long    0xffee42d1
  157.     .long    0x30273b76
  158. NaN:
  159.     .long    0x7fffffff
  160.     .long    0xffffffff
  161. p_Inf:
  162.     .long    0x7ff00000
  163.     .long    0x00000000
  164. m_Inf:
  165.     .long    0xfff00000
  166.     .long    0x00000000
  167.  
  168. ");    /* END ASM    */
  169. # endif    ERROR_CHECK
  170. __asm("
  171. .even
  172.     .globl _asin
  173. _asin:
  174.  
  175. ");    /* end asm    */
  176. #endif    /* __M68881__ || sfp004    */
  177. #ifdef    __M68881__
  178.  
  179.     __asm("
  180.     fasind    a7@(4), fp0    | asin
  181.     fmoved    fp0,a7@-    | push result
  182.     moveml    a7@+,d0-d1    | return_value
  183.     ");    /* end asm    */
  184.  
  185. #endif    __M68881__
  186. #ifdef    sfp004
  187.     __asm("
  188.     lea    0xfffa50,a0
  189.     movew    #0x540c,a0@(comm)    | specify function
  190.     cmpiw    #0x8900,a0@(resp)    | check
  191.     movel    a7@(4),a0@        | load arg_hi
  192.     movel    a7@(8),a0@        | load arg_low
  193.     movew    #0x7400,a0@(comm)    | result to d0
  194.     .long    0x0c688900, 0xfff067f8    | wait
  195.     movel    a0@,d0
  196.     movel    a0@,d1
  197.     ");    /* end asm    */
  198.  
  199. #endif    sfp004
  200. #if defined (__M68881__) || defined (sfp004)
  201. # ifdef    ERROR_CHECK
  202.     __asm("
  203.     lea    double_max,a0    |
  204.     swap    d0        | exponent into lower word
  205.     cmpw    a0@(16),d0    | == NaN ?
  206.     beq    error_nan    |
  207.     cmpw    a0@(24),d0    | == + Infinity ?
  208.     beq    error_plus    |
  209.     cmpw    a0@(32),d0    | == - Infinity ?
  210.     beq    error_minus    |
  211.     swap    d0        | result ok,
  212.     rts            | restore d0
  213. ");
  214. #ifndef    __MSHORT__
  215. __asm("
  216. error_minus:
  217.     swap    d0
  218.     moveml    d0-d1,a7@-
  219.     movel    #63,_errno    | errno = ERANGE
  220.     pea    _Overflow    | for printf
  221.     bra    error_exit    |
  222. error_plus:
  223.     swap    d0
  224.     moveml    d0-d1,a7@-
  225.     movel    #63,_errno    | NAN => errno = EDOM
  226.     pea    _Overflow    | for printf
  227.     bra    error_exit    |
  228. error_nan:
  229.     moveml    a0@(24),d0-d1    | result = +inf
  230.     moveml    d0-d1,a7@-
  231.     movel    #62,_errno    | NAN => errno = EDOM
  232.     pea    _Domain        | for printf
  233. ");
  234. #else    __MSHORT__
  235. __asm("
  236. error_minus:
  237.     swap    d0
  238.     moveml    d0-d1,a7@-
  239.     movew    #63,_errno    | errno = ERANGE
  240.     pea    _Overflow    | for printf
  241.     bra    error_exit    |
  242. error_plus:
  243.     swap    d0
  244.     moveml    d0-d1,a7@-
  245.     movew    #63,_errno    | NAN => errno = EDOM
  246.     pea    _Overflow    | for printf
  247.     bra    error_exit    |
  248. error_nan:
  249.     moveml    a0@(24),d0-d1    | result = +inf
  250.     moveml    d0-d1,a7@-
  251.     movew    #62,_errno    | NAN => errno = EDOM
  252.     pea    _Domain        | for printf
  253. ");
  254. #endif    __MSHORT__
  255. __asm("
  256. error_exit:
  257.     pea    _Error_String    |
  258.     pea    __iob+52    |
  259.     jbsr    _fprintf    |
  260.     addl    #12,a7        |
  261.     moveml    a7@+,d0-d1
  262.     rts
  263.     ");
  264. # else    ERROR_CHECK
  265. __asm("rts");
  266. # endif    ERROR_CHECK
  267.  
  268. #endif /* __M68881__ || sfp004    */
  269.